Cold neutrons trapped in external fields 
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The properties of inhomogeneous neutron matter are crucial to the physics of neutron-rich nuclei 
and the crust of neutron stars. Advances in computational techniques now allow us to accurately 
determine the binding energies and densities of many neutrons interacting via realistic microscopic 
interactions and confined in external fields. We perform calculations for different external fields and 
across several shells to place important constraints on inhomogeneous neutron matter, and hence 
the large isospin limit of the nuclear energy density functionals that are used to predict properties of 
heavy nuclei and neutron star crusts. We find important differences between microscopic calculations 
and current density functionals; in particular the isovector gradient terms are significantly more 
repulsive than in traditional models, and the spin-orbit and pairing forces are comparatively weaker. 
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The properties of inhomogeneous neutron-rich mat- 
ter are important in both astrophysical and terrestrial 
regimes. While the equation of state (EOS) and the 
pairing gap for homogeneous neutron matter have been 
studied extensively in microscopic theories [IHlj, inho- 
mogeneous neutron matter has received comparably lit- 
tle attention. Understanding the inner crust of neutron 
stars, which affects transient stellar cooling and deter- 
mine oscillation modes requires knowledge of inhomoge- 
neous neutron-rich matter [5H7]. Neutron-rich nuclei are 
also the subject of intense theoretical and experimental 
investigations, driven by their relevance for r-process nu- 
cleosynthesis as well as the intrinsic interest in the prop- 
erties of nuclei at large isospin; [5J [5] they are the princi- 
pal thrust of rare-isotope accelerators [ID]. 

Simulations of both the crust of neutron stars and of 
large neutron-rich nuclei employ nuclear energy density 
functionals fit to nuclei. These density functionals have 
proved to be extremely successful in describing many nu- 
clei, but involve large extrapolations to reach inhomoge- 
neous neutron matter. To test these extrapolations, we 
perform calculations of neutron drops - neutrons con- 
fined by artificial external fields and interacting via real- 
istic two- and three-nucleon forces. We vary substantially 
the number of neutrons as well as the strength and shape 
of the external fields to test the density functional. 

The EOS of homogeneous neutron matter has often 
been included as a constraint to density functional the- 
ories (eg. our objective is to allow inhomogeneous 
neutron matter to be employed in a similar manner. We 
find, for example, that once the bulk terms are fixed from 
the neutron matter EOS, the closed shells of neutrons are 
primarily sensitive to the gradient terms in the density 
functional. These pure neutron matter gradient terms 
have modest effects on nuclei, and hence they are not 
well constrained in fits to nuclear masses O H3] . The 
closed-shell systems are nearly independent of spin-orbit 
and pairing terms, but ground and excited states of a 
single neutron outside a closed shell, or of a single neu- 



tron hole, are a sensitive probe of the spin-orbit inter- 
action. Mid-shell results are sensitive to both spin-orbit 
and pairing terms. We compare our calculated results to 
several "standard" Skyrme models, and also to a model in 
which the isovector terms are adjusted to reproduce the 
ab-initio calculations; these changes are expected to have 
only a small effect on the nuclear energies used to fit the 
original parameters. The goal of these studies is to deter- 
mine which terms in the density functional can be probed 
through microscopic calculations, and how the adjusted 
values compare to traditional models. A realistic im- 
proved density functional will require a complete refitting 
of nuclear properties along with the properties of homo- 
geneous and inhomogeneous neutron matter [141 115] . 

Interaction and Methods: We report calculations of 
neutrons in harmonic oscillators (HO) of two frequencies 
and a Woods-Saxon (WS) well. The full Hamiltonian is: 
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where V t = (moj 2 /2)rf (HO) or V. t = -V /(l + exp[(r, - 
r )/a]) (WS) with V = -35.5 MeV, r = 3 fm and 
a = 1.1 fm, and h 2 /m = 41.44 MeV-fm 2 . The neutron- 
neutron potential is AV8' [IB], a slightly simplified 
version of the AV18 potential [17] : we find less than 0.25% 
differences in neutron-drop energies for these two poten- 
tials. We also add the Urbana IX model (UIX) [TBj three- 
nucleon interaction (TNI), including the p-wave two-pion 
exchange (Fujita-Miyazawa) TNI and a short-range phe- 
nomenological repulsion. We use this combination of two- 
and three-nucleon interactions because it produces an 
EOS consistent with known neutron star masses [T], and 
because several present-day Skyrme models have used 
this EOS to constrain the properties of homogeneous 
neutron matter. Further studies with different interac- 
tion models will be valuable, in particular to look at the 
spin-orbit interactions which might be increased with a 
three-pion exchange TNI as in Illinois- 7 |18j . 
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FIG. 1. (color online) Energies divided by hioN A ^ z for neu- 
trons in HO fields with hu = 10 MeV (top) and 5 MeV 
(bottom). Filled symbols indicate ab initio calculations; the 
dashed lines are Thomas- Fermi results (see text); the lower 
curves are from the SLY4 interaction and the upper curves 
show the modified SLY4 interaction described in the text. 
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FIG. 2. (color online) Energies per particle for neutrons in 
the Woods-Saxon field, symbols as in Fig. [T] 

Calculations are performed using Green's Function 
Monte Carlo (GFMC) [TS] and Auxiliary Field Diffu- 
sion Monte Carlo (AFDMC) [20] quantum Monte Carlo 
(QMC) methods. These algorithms evolve an initial trial 
state, ^t, in imaginary time to yield the ground-state. 
GFMC sums explicitly over spin and isospin states, and 
can use very sophisticated |16) . However it is lim- 
ited to small systems, up to 16 neutrons. In addition 
to sampling the spatial integrals as in GFMC, AFDMC 
also samples the spin and isospin degrees of freedom, and 
hence it can treat larger systems [3]. Both methods use 
a constraint involving the overlap with to eliminate 
the Fermion sign problem, and hence are approximate. 
Studies of light nuclei and neutron matter show they give 
results within 1% of the exact ground-state energy. 

We use external fields yielding low or moderate den- 
sities. However, even at small densities neutrons are 
strongly interacting and pairing can be important. Re- 
cent microscopic calculations of neutron matter give s- 
wave pairing gaps of several MeV [H [ST] . One- and two- 
nucleon properties including pairing gaps and spin-orbit 
splittings can be more sensitive to models of the three- 
nucleon interaction. Calculations of very small neutron 
drops (N=6,7,8) have been performed previously [22H24] . 
Even these calculations indicated a substantial differ- 
ence with traditional Skyrme models, which overbind the 
drops and give too-large spin-orbit splitting. 
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TABLE I. Selected energies from GFMC and AFDMC cal- 
culations using AV8'+UIX with HO external fields. 



Results: The ground-state energies versus neutron 
number N for the HO potentials are given in Fig. [T] and 
for the WS potential in Fig. [2] Up to N=16 both GFMC 
and AFDMC results are included. They agree very well 
for the 10-MeV HO interaction, while for ftio — 5 MeV, 
the AFDMC results are slightly higher than the GFMC; 
the maximum difference is 3%, and more typically results 
are within 1%. The bigger difference for the lower density 
5-MeV drops presumably arises because the AFDMC 
does not yet include pairing, while the GFMC does. 

In addition to the microscopic calculations, re- 
sults for several different Skyrme models are shown 
in Fig. [T] We also show results for Thomas-Fermi 
local density approximations [25] using E(p n )/N — 
^(3/5)(fi, 2 /2m)(37r 2 p n ) 2 / 3 ; the upper horizontal line is for 
free particles, £ = 1, and the lower has £ = 0.5, a reason- 
able approximation to the EOS of low-density neutron 
matter. For the 10-MeV well, the density functionals 
give energies significantly below the Monte Carlo results 
for all N. The energies are also lower for the 5-MeV well, 
but less so. This overbinding is a general feature of all 
the Skyrme models considered. It is intriguing that these 
same Skyrme models underbind the properties of very di- 
lute neutron systems, typically they are fit to the neutron 
matter EOS at p = 0.04 fm -3 and above. 

Since the Skyrme homogeneous neutron matter EOS 
have been fit to various microscopic calculations, this 
overbinding suggests that the gradient terms in inhomo- 
geneous neutron matter should be more repulsive. The 
observed differences between ab-initio results and the 
Skyrme functionals are much larger than the differences 
between experiments and Skyrme models in nuclei, as 
expected, because of the large extrapolations to inhomo- 
geneous neutron matter. 

Isovector Gradient Contributions: As is apparent in 
Fig. [T] for harmonic oscillators there are closed shells 
at N= 8, 20, and 40 neutrons. These closed-shell states 
are almost exclusively sensitive to the neutron matter 
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EOS and the isovector gradient terms; pairing and spin- 
orbit play nearly no role. Hence they are direct probes of 
the gradient terms; to examine them we have altered the 
isovector gradient terms in the SLY4 interaction [TT] to 
approximately reproduce the QMC results using a mod- 
ified version of the ev8 code , The gradient terms are 
adjusted without changing any isoscalar (T=0) parame- 
ters or the homogeneous neutron matter EOS. 

The lowest-order gradient contribution to the energy 
density for inhomogeneous matter is Gd[Vp n } 2 . The con- 
stants Gd are small and often negative, for example, 
G d = -16, -7, 17, -17, -7 MeV-fm 5 for the SLY4, SLY7, 
BSK17, SkM*, and SkP interactions. Repulsive gradient 
terms for neutron matter are to be expected on rather 
general grounds, and are required for the absolute sta- 
bility of uniform matter in the absence of a background 
field. The adjusted interaction SLY4-adj gives Gd = 26.5, 
a similar adjustment to the BSK17 interaction which is 
more attractive for homogeneous neutron matter yields 
G d = 64. A single adjustment of Gd markedly improves 
the agreement with QMC results for both the HO and 
WS fields. A precise fit to both neutron matter and these 
results would require a more general form of the density 
functional. 

Isovector Spin- Orbit: By examining neutron numbers 
slightly away closed shells, we can constrain the spin- 
orbit interaction for neutron drops. For example, N — 
7, 9 results are sensitive to the spin-orbit interaction, but 
not to the pairing terms. We find the spin-orbit interac- 
tion to be small in the calculated drops; the energies for 
some low-lying spin-orbit partners are given in Table [I] 

Small spin-orbit splitting had been found previously in 
calculations of N=7 drops and our results show similar 
effects for all systems near closed shells (N=7,9,ll,...). 
The simplest (standard) Skyrme parametrizations give a 
strength ratio of 3:1 between isoscalar (T = 0) and isovec- 
tor (T = 1) spin-orbit couplings. We find an even smaller 
isovector coupling, approximately 1/6 of the isoscalar 
coupling, reproduces the ab initio calculations. The com- 
bined factor of 1/6 is in reasonable agreement with 1/8 
found in a diagrammatic examination of spin-orbit split- 
tings from microscopic interactions |27j . and with re- 
sults obtained in an examination of Skyrme parameters 
from Pb isotopes [55]. Relativistic mean fields yield zero 
strength in the isovector channel [29 , while recent results 
from the pion contribution from chiral interactions give 
nearly equal isovector and isoscalar strengths [30 . 

Isovector Pairing: The mid-shell results (eg. N=14, 
30) and odd-even staggerings are sensitive to the pairing 
interactions as well as the spin-orbit force. Fixing the 
spin-orbit strength from near closed-shells, we adjust the 
pairing strength to fit the calculated spectra. There is a 
significant interplay between the pairing and spin-orbit 
forces required to reproduce microscopic calculations. A 
small spin-orbit force results in many quasi-degenerate 
levels which enhances pairing in mid-shell systems. 
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FIG. 3. (color online) Calculated radii of neutrons confined 
in HO (upper) and WS (lower), fields compared to original 
and adjusted Skyrme models (see text). 

Several models of pairing are used in density- functional 
theories. We employ a simple volume parametrization 
with a delta-function spatial dependence, a density cut- 
off that restricts pairing to p n < Pa, and limit the pairing 
to single-particle orbitals less than 5 MeV from the Fermi 
energy. We find a reduction from a typical 1 GeV-fm 3 
strength to half that value significantly improves agree- 
ment with microscopic results. A reduction of pairing 
in neutron-rich nuclei has recently been found to give a 
better fit to experimental energy differences of 156 nuclei 
of mass A=118 to 196 [51] . 

Adjusting these three parameters (gradient term Gd 
= 26.5, spin-orbit coupling = 123 MeV-fm 3 and pairing 
strength = 500 MeV) in the density functional increases 
the agreement across all external fields and all particle 
numbers. This is shown by the upper solid curves (SLY4- 
adj) in Fig.[l] 

Radii and Mass Distributions: Our calculations yield 
precise estimates for RMS radii and the density distribu- 
tions of the smaller drops. The average densities, defined 
as J d 3 rp 2 l (r)/N, of the drops in the 5-MeV HO well are 
approximately 0.02 fm~ 3 , or about 1/8 nuclear matter 
saturation density, while for the 10-MeV HO and the 
WS wells they are ~ 0.045 fm~ 3 , or almost 1/3 nuclear 
matter saturation density. 

The RMS radii obtained in microscopic calculations 
are compared with the original and adjusted Skyrme den- 
sity functional results in Fig. [3] The density distribu- 
tions for N=8 and 14 are compared in Fig. [4] Since we 
are comparing gross properties of inhomogeneous mat- 
ter, we plot the densities weighted with the phase space: 
r 2 p n (r), which gives a better picture of the density distri- 
butions near the average density of the system. In every 
case the adjusted Skyrme interaction produces a better 
description of the radii and density distributions. The 
N = 8 calculations depend primarily upon the gradient 
terms, the reduction in pairing and spin-orbit are also im- 
portant for N = 14. The improvement in the mid-shell 
N=14 case is particularly dramatic, as a significant shift 
in the density occurs with the modified isovector Skyrme 



4 



0.3 

°- 0.1 
u 



,--—0.4 

& 03 
clO-2 
r V 0.1 





N=8 

/ ^^^^^^ • - _ 


■ to = 5 MeV 

■ to = 10 MeV 

SLY4-adj 

-- SLY4 




i , i , i , i , T 







3 

r (fin) 



FIG. 4. (color online) Calculated densities of neutrons in HO 
potentials, compared to Skyrme models (see text). 



parameters, bringing the results into much better agree- 
ment with microscopic calculations. 

Conclusions: We have examined the properties of neu- 
trons confined in external fields to study the properties 
of inhomogeneous neutron matter. These ab-initio calcu- 
lations place significant constraints on the nuclear energy 
density functional in a regime far from that probed by fit- 
ting to available nuclei. They indicate the need for more 
repulsive gradient terms in pure neutron matter, and a re- 
duced isovector spin-orbit and pairing strength compared 
to standard functionals. With a combined fit of density 
functionals to both nuclei and neutron matter, more reli- 
able predictions should be possible for very neutron-rich 
nuclei including those participating in r-process nucle- 
osynthesis. These improved functionals would also be ex- 
tremely valuable in examining the neutron skin thickness 
of lead 32J , as can be probed in parity- violating electron 
scattering. Much more reliable predictions for extremely 
neutron-rich astrophysical environments can also be ex- 
pected. The numerical values of the results shown in the 
figures are given in Ref. |33j . 
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